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ABSTRACT 

Let xi , • • ■ , xt G K" . A simultaneous integer relation (SIR) 
for Xi, ■ • ■ , Xt is a vector m £ Z" \ {0} such that Xi"^m — 
for i — I,-- - ,t. In this paper, we propose an algorithm 
SIRD to detect an SIR for real vectors, which constructs an 
SIR within 0{n^+n'^ log A(X)) arithmetic operations, where 
X{X) is the least Euclidean norm of SIRs for xi,--- ,xt. 
One can easily generalize SIRD to complex number field. 
Experimental results show that SIRD is practical and bet- 
ter than another detecting algorithm in the literature. In 
its application, we present a new algorithm for finding the 
minimal polynomial of an arbitrary complex algebraic num- 
ber from its an approximation, which is not based on LLL. 
We also provide a sufficient condition on the precision of the 
approximate value, which depends only on the height and 
the degree of the algebraic number. 

1. INTRODUCTION 

Let xi, ■ • ■ , Xt be vectors in R", and denote (xi, ■ • ■ , xt) 
by X. A simultaneous integer relation (SIR) for xi, • ■ • ,xt 
is a vector m G \ {0} such that X^m = 0, i.e. Xi"^m — 
for i — 1, • • • ,t. For short, we also call m an SIR for X. 
When t = 1, we say that m is an integer relation for xi. 
The problem of detecting integer relations for a rational or 
real vector is quite old. Historical surveys can be found in 
O 1141 1101 1171 113j . Among these integer relation detecting 
algorithms, the HJLS algorithm |16l 117) and the PSLQ al- 
gorithm [121 113] have been used frequently. 

In the present paper, using the technique to construct the 
hyperplane matrix in HJLS and a generalized method of 
the matrix reduction from PSLQ we propose an algorithm 
SIRD, which can be used to detect an SIR for t real vectors. 
The cost of our algorithm is at most ©(n* -I- log A(X)) ex- 
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act arithmetic operations for detecting an SIR for X, where 
X{X) represents the least Euclidean norm of SIRs for X. 
Furthermore, our detecting algorithm SIRD either always 
finds an SIR for X if one exists or proves that there are no 
SIRs for X of norm less than a given size. Experimental 
results show that SIRD is practical. 

In application, we successfully apply SIRD to find the 
minimal polynomial of an algebraic number a G C with 
degree and height at most n and H respectively from its 
an approximation a satisfying maxi<i<„ ]a* — a'| < e, and 
propose the corresponding algorithm MPF, where the min- 
imal polynomial of an algebraic number a is the unique 
primitive polynomial p{x) G Z[x] of least degree such that 
p{a) = 0. In fact, for i from 1 to n we run SIRD with vi = 
(l,Re(Q),--- ,Re(a^))^, V2 = (0, Im(a), ■ • ■ , Im(a*))^ as 
its input and then an exact SIR for vi , V2 has been de- 
tected. We provide a sufficient controlling on e and prove 
that such an e is sufficient to enable an exact SIR for vi 
and V2 to be also an SIR for (l,Re(Q),--- ,Re(a'))"^ and 
(0, Im(a), ■ • ■ , Im(a'))"^, where e depends only on n and H, 
as in (|5.3|l . It implies the correctness of MPF and is better 
than already existing results in [181 [25] . 

1.1 Related Works 

In 16l Il7j. J. Hastad, B. Just, J. C. Lagarias, and C. P. 
Schnorr not only presented the HJLS algorithm and the first 
rigorous proof of a 'polynomial time' bound for a relation 
finding algorithm but also proposed a simultaneous relations 
algorithm (see [171 section 5]), whereas HJLS is numerically 
unstable. The unstable examples can be found in [121 113| . 
In their draft [26[ . C. Rossner and C. P. Schnorr studied the 
case of t = 2 by using a modified HJLS algorithm. But for 
the moment, [26] is still in a preliminary state with some 
open problems. The PSLQ algorithm, together with related 
lattice reduction schemes such as LLL [21], was named one 
of ten "algorithms of the twentieth century" by the publi- 
cation Computing m Science and Engineering (see [111 1^). 
and is now extensively used in Experimental Mathematics, 
with applications such as identification of multiple zeta con- 
stants, a new formula for tt, finding algebraic relations and 
so on (see [4l|3][2]). Moreover, PSLQ is numerically stable 
and can be easily generalized to complex number field and 
Hamiltonian quaternion number field (see [l3]), but it is not 
suitable to detect an SIR for several real vectors. 



The SIRD algorithm in this paper is to detect an SIR for t 
real vectors and can be applied to detect an integer relation 
in for a complex vector or a Hamilton quaternion number 
vector. A significant body of experimental data shows that 
SIRD is practical and better than the HJLS simultaneous 
relations algorithm. 

In fact, the MPF algorithm in this paper is a positive 
answer to the following interesting question: Suppose we 
are given an approximation to an algebraic number a, and 
two bounds on the degree and the size of the coefficients 
of its minimal polynomial respectively. Is it possible to in- 
fer the minimal polynomial? The question was raised, in- 
dependently, by Manuel Blum in theoretical cryptography 
(see |19l [20] ) and the last author of this paper in automated 
reasoning (see [H]). The first complete answer to this ques- 
tion, KLL algorithm, was presented by R. Kannan, A.K. 
Lenstra and L. Lovasz in [191 120) by using the celebrated 
lattice reduction algorithm LLL [21]. In the computer alge- 
bra system Maple, the built-in function PolynomialTools : - 
MinimalPolynomialO is a function to find a polynomial of 
degree n (or less) with small integer coefficients which has 
the given approximation r of an algebraic number as one of 
its roots and is based on KLL algorithm. The correctness 
of the polynomial returned by the built-in function depends 
on the accuracy of the approximation (see Maple's Help). 
From another aspect, the minimal polynomial of an alge- 
braic number a with exact degree n can be found by detect- 
ing an integer relation for the vector v = (l,a, ••• ,a")'^. 
Besides HJLS, B. Just also presented an algorithm to detect 
integer relations for a given vector consists of algebraic num- 
bers in We can apply Just's algorithm or HJLS to the 
vector V for finding the minimal polynomial of a. However, 
both Just's algorithm and HJLS are not numerically stable, 
as mentioned previously. All these algorithms are based on 
LLL. Two authors of this paper presented a method to re- 
construct a rational number from its an approximation by 
using continued fraction in [30]. It may be viewed as an 
answer to a special case of the question. Based on PSLQ, 
one can find algebraic relations, such as [BlIZllSKI]) whereas 
these articles did not involve the minimal polynomial find- 
ing. The authors of this paper also presented an algorithm 
in [25] for finding the minimal polynomial of a real algebraic 
number from its an approximation. However, these PSLQ 
based algorithms can not deal with complex algebraic num- 
bers since PSLQ only outputs a relation in Gaussian integer 
ring for a complex vector. 

Fortunately, our simultaneous integer relation detection 
algorithm SIRD in present paper can be used to overcome 
these pitfalls. Applying SIRD to one or two real vectors, we 
present another affirmative answer, the MPF algorithm, to 
the question above. We show that MPF is a more efficient 
minimal polynomial finding algorithm comparing with the 
algorithms in [181 125[ and provide a sufficient condition on 
the error controlling, from which we can claim that the poly- 
nomial returned by MPF is the exact minimal polynomial 
of the algebraic number that we only know an approximate 
value and two bounds on its degree and height. Although 
a similar even better complexity can be obtained by KLL, 
MPF has its own meaning since it is a new method without 
using LLL reduction. 

Road-map. In section[2]and[2]we first give some prelimi- 
naries, and then present the SIRD algorithm and analyze it. 
We report on some experimental results about the perfor- 



mance of SIRD in section [IJ apply SIRD to find the minimal 
polynomial of an algebraic number from its an approxima- 
tion and propose the MPF algorithm in section [S] in which 
we also analyze MPF and present the result of error control- 
ling. We conclude this paper with section [S] 

Notations. Throughout this paper, Z, R, and C are the 
sets of integers, real numbers, and complex numbers respec- 
tively. The real and imaginary parts of z G C will be denoted 
Re{z) and Im(2:) respectively. For c G R, [c] = [c |J . 
All vectors in this paper are column vectors, and will be 
denoted in bold. If x £ R", then |lx||2 represents its Eu- 
clidean norm, i.e. |jx|j2 = y< x, x >, where < *, * > is 
the inner product of two vectors. We denote n x n iden- 



tity matrix by /„. Given a matrix A = 



), we denote 



its transpose by , its trace by tr{A), its determinant by 
\A\, and its Frobenius norm by \\A\\f = {ti{A'^ A))^^'^ , i.e. 
11^1!^' ~ {'12^1,])^^'^ ■ We say that a matrix A is lower trape- 
zoidal if J = for i < j. GL(n, Z) is the group of n x n 
unimodular matrix with entries in Z. The height of a vector 
is defined by the maximum of all the absolute values of its 
entries. For a polynomial f(x) = "^"^ofiX^, we denote by 
deg(/) its degree with respect to x, \\f\\i — X^r^o 1-^*1 
one norm, 

II/II2 = (Er=o \MY^^ its Euclidean length, and 
height(/) = maxo<i<n its height. 

2. PRELIMINARIES 

In what follows we always suppose that xi , • ■ • , xt are lin- 
early independent vectors in R", where Xi = {xi^i, • • • , Xi^n)"^ . 
Obviously, we have t < n. We denote by X the matrix 
(xi, • • • ,xt), and suppose that X £ R"'^* satisfies 



Kl,„_t + 1 X2,n-t + l 
Xl.n-t + 2 X2,n-t + 2 



Xl- 



X2,- 



Xt,n-t+\ 
Xt,n-t + 2 



Xt,n 



(2.1) 



unless otherwise specified. For X G R"^' not satisfying 
(|2.ip . exchanging some rows of X produces X' = CX, where 
C is an appropriate matrix in GL{n, Z). And then we detect 
an SIR for X' . If m is an SIR for X' , then C'^m is an SIR 
for X. 

2.1 Hyperplane Matrix 

Definition 2.1 (Hyperplane Matrix). Let X = {jci, 
■ ■ ■ ,xt) G R"*^*. A hyperplane matrix with respect to X is 
any matrix W G R"x("-') such that X'^W = and the 
columns ofW span X^ = {y G R" : Xi^y = 0,i = 1, ■ • • 

Now we introduce a method to construct a hyperplane 
matrix for X. 

Let bi, • • • , bn form a standard basis of R", i.e. the i-th 
entry of bi is 1 and others are 0. By performing the process 
of standard Gram-Schmidt orthogonalization to xi, • ■ • , xt, 
bi , ■ ■ • , bn in turn we have 



Xl 



Xk 



bl' 



Xl 

I|X1||2^ 

Xk* 

||Xk*||2- 



E<Xk,Xj > , 
Xj , 
, = 1 <Xj*,Xj* > 

1 ^ ~ 2, • ■ • , t. 



, ■r-^ < bi,Xj* > * , * bi 



< Xi*,Xi* > 



bi* 



2. Pjf 



> 



where /Xi,. 



bi* 



I|bi*|l2^ 

r < bi,bj* > 

<bj*,bj* >' 

lo, 



2, •■• ,71, 



if libi 



if l|bj 



= 0. 



Lemma 2.2. iet Xk, Xk*, bj and bj* 6e as above. Then 
1. there exist t elements tn{l, - ■ ■ ,n} denoted by ji - ■ ■ , jt 



such that bjj^ " 

2. bn-t + l* 



bj,* 



= 0. 
0. 



Proof. Part 1 easily follows from the process of standard 
Gram-Schmidt orthogonalization. We next 
• • • = bn* = when pTTj) holds. Set 



aixi H h atXt + hhi 



+ In-thn-t = 0. 



Taking each side as a column vector and observing the last 
t components of two sides, we have ai — ■ ■ ■ — at — 0. And 
since bi,--- , bn-t are linearly independent, we have h = 
• ■ ■ = In-t = 0. Thus the n vectors xi , ■ • ■ , xt , bi , ■ • ■ , bn-t 
are linearly independent. This implies that bn-t+i* = • ■ • = 

bn* =0. □ 



Definition 2.3 (Hx)- For X e K"""* satisfying h2l]) . 
define Hx to he the n x {n ~ t) matrix (bi*, • ■ • , bn-t*). 



Lemma 2.4. Let X e 

1. HxHx = In-t- 



and Hx be as above. Then 



x^ 

2. ||-ffx||F = y/n - t. 

3. (xi*, • • • , Xt*, Hx) is an orthogonal matrix. 

4. X^Hx =0, I.e. Hx is a hyperplane matrix of X . 

5. Hx IS a lower trapezoidal matrix and every diagonal 
element of Hx is nonzero. 

Proof. Since every two columns of Hx are orthogonal, 
part 1 follows. And part 2 follows from part 1. Let X* — 
(xi*,--- ,xt*)^. Obviously, (xi*,--- ,xt*,Hx) is an or- 
thogonal matrix. From part 3 and standard Gram-Schmidt 
orthogonalization we have X*^ Hx = and X — X*Q re- 
spectively, where Q is an appropriate t x t invertible ma- 
trix. Thus X^ Hx ~ X*^ Hx ~ and hence that part 
4 follows. We now prove part 5. Denote the fc-th ele- 
ment of bi* by 6* J.. The diagonal elements of Hx are 
b*i for i — 1, - ■ ■ ,n — t. Before normalizing bi* we have 
6* J = 1 — X]fc=i ^fc^i ~ X]j=i^j,^!i ^iid at the same time, 

/ ||bi* Hi bi*, bi* 1 - ELi - Ei^l Thus 

all the diagonal elements of Hx are nonzero. Now we only 
need to show that Hx is lower trapezoidal. From stan- 
dard Gram-Schmidt orthogonalization, we can check that 
bi.k =< bi*,bk* holds for i > k. This completes the 
proof. □ 

So far, we have had a method to produce a hyperplane 
matrix ii"x for X G R"""*. The basic idea is from HJLS 
(see [161 117j). The same strategy was also used in PSLQ, 
however, in which partial sum was adopted instead of Gram- 
Schmidt orthogonalization. 



Lemma 2.5. For X 
HxH]^. Then 
1. P^ = Px. 



(xi. 



,xt G 



define Px 



In X^i — 1 Xi Xi 



I \\Px\\f = v^fr^. 

5. Pxz — z for any z G X''^ . Particularly, Pxm — m for 
any SIR m for X . 

Proof. The proof of the first part is easy. Let U = 
(xi*, • • • , Xt*, Hx). From Lemma [2.41 we have /„ = UU^ = 
HxHx + X^i^i Xi'xi*"^. Thus part 2 follows. Part 3 and 
part 4 follow' from P^ = Hx{H'^Hx)H'^ = HxH'^ = Px 
and ||Px|||. = tr(Pl'Px) = tr(Px) = tr(//^//x) = n~t 
respectively. Since z G X^ , we have < Xi,z >= for 
i = 1- ■ ■ ,t. And the process of standard Gram-Schmidt 
orthogonalization implies Xi*"^z — 0. Thus we have Pxz = 
z — Xi*Xi*"^)z — z from part 2. □ 

From Lemma [2.4l and Lemma [2. 5 1 we can easily generalize 
the Theorem 1 in [13] to the case of X G R"""*. 

Theorem 2.6. Let X G R"^* and Hx be as above. Sup- 
pose that for any matrix A G GL{n,'Z) there exists an or- 
thogonal matrix Q G r("-*' ><("-*) .such that (hij) = AHxQ 
is lower trapezoidal and all of the diagonal elements of {hij) 
satisfy hjj 7^ 0. Then for any SIR m of X we have 

' ' • ' (2.2) 



maxi<j<„-t \hj 



mm -r— 

l<j<n-t \hj 



< m 



As this theorem easily follows from the proof of Theorem 1 
of [13J with little modifications, the detail has been omitted 
here. 

The lower bound given in ([2.2p when t = 1 is consistent 
with a similar lower bound in [14II15| . Moreover, if a method 
to reduce the norm of Hx by multiplication by some uni- 
modular A G GL{n, Z) on the left has been developed, then 
it will produce an increasing lower bound on \(X), where 
X{X) is the least Euclidean norm of SIRs for X. In fact this 
theorem suggests a strategy to detect an SIR for X. 

2.2 Matrix Reduction 

We now study how to reduce the hyperplane matrix Hx- 
First we recall (modified) Hermite reduction in [13[ . 

Algorithm 1 (Modified Hermite Reduction). 

Input: a lower trapezoidal matrix H = [hij) G R"*^'""^) 

with hjj 7^ 0. 
Output: a reducing matrix D of H. 
1: D := In 

2: for i from 2 to n do 

3: for j from i — 1 by —1 to 1 do 

4: q := [hij/hj^j], where [c] — [c+ 1/2J for a real 

number c. 
5: for k from 1 to n do 

6: di.fc :~ di^k — gdj,fc 

7: return the n x n matrix D. 

If Algorithm [T] output D for an n x (n — 1) matrix H, 
we say that DH is the modified Hermite reduction of H and 
that D is the reducing matrix of H . This reduction develops 
the left multiplying modified Hermite reducing matrix D. 

Hermite reduction is also presented in [13], and is equiv- 
alent to modified Hermite reduction for a lower triangular 
matrix H with hjj 7^ (see [131 Lemma 3]). Both the two 
equivalent reductions have the following properties: 



1. The reducing matrix D € GL{n,1,). 

2. For all A; > i, the (modified) Hermite reduced matrix 
H' = (M,,) = DH satisfies < \hi^,\/2 = \h.u\/2. 

In order that the reduced and reducing matrices of Hx G 
]]jnx(n-t) gg^^jgfy ^jjg properties above, we need the fol- 
lowing generalized Hermite reduction. 

Algorithm 2 (Generalized Hermite Reduction). 



Algorithm 3 (The SIRD Algorithm). 



Input: a lower trapezoidal matrix H = (hij) £ 

with hjj 7^ 0. 
Output: a reducing matrix D of H . 

1 

2 

3 



X (ri — t) 



1 



9: 
10: 



D ~I„ 

for i from 2 to ti do 

if i < n — t + 1 then temp := i - 
else temp := n — t 
for J from temp by —1 to 1 do 
Q ■= [hij/hjj] 
for k from 1 to n do 

for every two integers si, S2 G {n — t + 1, ■ ■ ■ ,n} satis- 
fying si < S2, /isi,n-t = and hs2.n-t / do 

exchange the si-th row and the S2-th row of D. 
return the n x n matrix D. 



If Algorithm [2] output D for an n x {n — t) matrix H , 
we call DH the generalized Hermite reduction of H and _D 
the reducing matrix of H . Obviously, generalized Hermite 
reduction is equivalent to modified Hermite reduction when 
t = 1. In addition, we can easily check that generalized Her- 
mite reduction remains the two properties mentioned above. 



Input: (xi, • • • ,xt) = A G R"^* satisfying (gH]) 
Output: either output an SIR for X or give a lower bound 
on \{X). 

1: Initiation. Compute the hyperplane matrix Hx, set 
H := Hx, B := In- 

2: Reduction. Call Algorithm [2] to reduce Hx produc- 
ing the reducing matrix D G GL{n,Z). Set X : = 
XD-^,H := DH, B := BD'^. 

3: loop 

4: Exchange. Let H — (hi.j). Choose an integer r such 
that 7''|/ir_^| > 7*|/ii,i| for 1 < i < n — t, where 7 > 
2/%/3. Define the permutation matrix R to be the 
identity matrix with the r and r + 1 rows exchanged. 
Update X := XR, H ~ RH, B := BR. 

5: Comer. Let 



a :— hr,r, 

X :— hr + l,r + l. 



P :— hr + l,r, 

5:= TFT 



A2. 



(3.1) 



and r + 1 be 



Update H ~ HQ. 



Let Q := In-t. If r < 71 — t, then let the submatrix of 
Q consisting of the r-th and (r + l)-th rows of columns 

fp/s -A/^^ 
{x/s p/s 

6: Reduction. Call Algorithm [2] to reduce Hx producing 
D. Update X := XD'^H ~ DH, B ~ BD~\ 

7: Compute G :— 1/ maxi<j<„_t|/ij,j|. Then there ex- 
ists no SIR whose Euclidean norm is less than G. 

8: if Xj = for some 1 < j < n, or h„-t,n~t = then 

9: return the corresponding SIR for X. 

10: end loop 



Remark 1. There are two main differences between (mod- 
ified) Hermite reduction and generalized Hermite reduction. 
Firstly, the last t — 1 rows of H will also be reduced by 
the first n — t rows of H in generalized Hermite reduction, 
while (modified) Hermite reduction can not do so. Sec- 
ondly, generalized Hermite reduction exchanges the si-th 
row and the S2-th row of D if si < S2, hs^,n-t ~ and 
hs2,n-t 7^ (from Step |5] to Step This implies that 
if h„-t+i,n-t ~ after generalized Hermite reduction then 
h„-t+2,n-t = • • ■ = h„,n-t = 0. This property plays an 
important role in the proof of Lemma |3. II 

3. THE SIRD ALGORITHM 

3.1 The Description of SIRD 

Using the hyperplane matrix constructing method and 
generalize Hermite reduction in the previous section we can 
get a simultaneous integer relation detecting algorithm SIRD. 

3.2 Analysis of SIRD 

Let H{k) be the result after k iterations of SIRD. 

Why do we set the parameter 7 > 2/y/3 at Step|4]? Sup- 
pose the r chosen in Step |4] is not n ~ t. In this case we let 
a, /3, A, 5 be as in (|3T|) . Then 

a 
/? A 

is the submatrix of H{k — 1) consisting of the r and r + 1 
rows of columns r and r -\- 1, where r < n — t. After Step 
|4]has been performed A may not be zero, which makes that 



H is not lower trapezoidal. After Step [5] the result is 



fP/5 -X/5 
\X/5 p/5 



S 

aj3/5 





-aX/5 



(3.2) 



Since r is chosen such that 7'"|/ir,r(fc — 1)| is as large as 
possible, and r < n— twehave |ftr+i,r+i(fc— 1)| < ^\hr,r{k~ 
1)1, hence |A| < ^\a\. From the property of generalized 
Hermite reduction we have that ]/3] < which then gives 

5 liWTxF 



r,r (fc — 1) 



< 



1 1 



T 
< 1, i.e. 



(3.3) 



Thus \hr,r\ is reduced as long as + < 1, i-C- 7 > 

2/\/3. As was pointed out by Borwein (see [8]), although 
this increases this is not a significant problem. At 

each step we force the larger diagonal elements of H toward 
h„^t,n-t, where their size can be reduced by at least a factor 
of 2 when r — n — t. 

As a matter of fact, the parameter 7 can be freely chosen 
in the open interval (2/\/3, +00). 

Lemma 3.1. // hjj{k) — for some 1 < j < n — t and 
no smaller k, then j — n — t and an SIR for X must appear 
as a column of the matrix B. 

Proof. By the hypothesis on k we know that all diagonal 
elements of H{k — 1) are not zero. Now, suppose the r 
chosen in Step |4] is not n — t. Since generalized Hermite 
reduction does not introduce any new zeros on the diagonal, 
and from the analysis of Step [4] and Step [5] above, we have 



that no diagonal element of H{k) is zero. This contradicts 
the hypothesis on k and our assumption that r < n — t was 
false. Thus we have r = n — t after the (fc — l)-th iteration 
has been completed. 

Next we show that there must be an SIR for X appeared as 
a column of the matrix B. We have Hx = from Lemma 
[lHand hence that = X'^BB-^Hx = X'^BB-^HxQ = 
X'^BH{k~^ 1), where Q is an appropriate orthogonal (n — 



i) X (n — t) matrix. Let (zi, ■ • ■ , zt) 
• ■ ■ , Zi^nY" . Then 



— X B, where Zi = 



Vo 



= X^BH{k- 1) = 




H{k- 1) 



Zl,n~th„~t,n~t{k — 1)\ 
Zt,n — thn-t,n-t{k — l)y 

We know /i„_t+i,n-t(fc— 1) = and /i„_t,„_f (fc— 1) / from 
hn-t,n-t{k) = 0. From Remark[T]and /i„-t+i,„-t(fc— 1) = 
we have hn-t+2,n-t{k ~ I) = ■ ■ ■ = hn,n-t{k — 1) = which 
implies the last equality. Since hn-t,n-t{k — l) 7^ 0, it follows 
that Z\,n-t — ■ ■ ■ = zt,n~t ~ 0. Thus the {n — t)-th column 
of B is an SIR for X.' □ 



From Theorem 12.61 and Lemma 13.11 the correctness of 
SIRD has been proved. Moreover, we have 

Theorem 3.2. Let \(X) be the least Euclidean norm of 
any SIR for X . Let m be an SIR detected by SIRD. Then 
l|m||2 < 7""'~U(X) for all 7 > 

Proof. Assume r — n — t with h„-t,n-t{k) 7^ and 
hn-t,n-t(k + 1) = at the fc-th iteration of SIRD. Then 
from Theorem 12.61 and the exchange rule of SIRD we have 



A(X) > 1/ max i/ii,i(fc)i > 7 

Ki<n — t 



t{k)\. 



At this time, ||m||2 = l/\hn-t,n-t{k)\ holds from the same 
strategy in the proof of Lemma 10 in [13] . □ 

Definition 3.3 (the II function). For the k-th iter- 
ation in SIRD, define 



n(fc) = Y\_ ™n 

l<j<n-t 



7 



'A(X), 



\hjAk)\ 



The routine of analyzing the number of iterations in [13] 
can be carried over here with redefining the 11 function as 
above. So we state the following lemma directly without 
proof. 



Lemma 3.4. For k > 1 we have 



1. (7"-'A(X)) 
the least norm of SIRs for X . 

2. n(fe)> y^n(fc-i). 



> n(fc) > 1, where \{X) is 



From this lemma, it follows that the II function is increas- 
ing with respect to k and has an upper bound for a fixed 
7 £ (2/\/3, -l-cxa). Thus we have 



Theorem 3.5. If X £ R"'^* has SIRs, then the number 
of iterations such that SIRD finds an SIR for X will be no 
more than 



log(7"-*A(X)) 



Proof. From Definition 13.31 we can infer 11(0) > 1. And 
by Lemma 13.41 we know that 



U{k) > 



472 



72+4 



Solving k from this inequality gives the conclusion, as was 
to be shown. □ 

Corollary 3.6. IfX e R"^' has SIRs, then there exists 
a 7 such that SIRD will find an SIR for X in polynomial 
time 0{n'^ + n^ logA(X)). 

Proof. Let 7 = 2. Then SIRD wiU construct an SIR for 
X in no more than 

(n - tf{n + t-l) + {n~- t){n + t- 1) log A(X) 

iterations. SIRD takes 0{n — t) exact arithmetic operations 
per iteration, and hence that 0((n — t)^ + (n — t)^ log \{X)) 
exact arithmetic operations is enough to produce an SIR for 
X. Since t < n, the proof is complete. □ 

Remark 2. From this corollary, we can claim that our de- 
tecting algorithm always return an SIR for X if one exists. 
Additionally, SIRD will produce lower bound on the Eu- 
clidean norm of any possible SIRs for X (Theorem 12. 6p . 
Thus SIRD can be used to prove that there are no SIRs for 
X of norm less than a given size. 

Remark 3. PSLQ may be viewed as a particular case of 
SIRD when t = 1. Similarly with PSLQ, SIRD can be eas- 
ily generalized to complex field with 7 > \/2 such that the 
outputs are in Gaussian integer ring and all conclusions men- 
tioned above hold with corresponding modifications. 

Remark 4- Moreover, SIRD can also be applied to detect 
an integer relation in for a given complex vector. For 
example, suppose z = x4-y7 in C" with vector components 
X, y e R" where I = Then SIRD can give an SIR m 

for (x,y), and hence that m £ is an integer relation for 
z, but PSLQ only can give a Gaussian integer relation in 
Z[/]". This is one of the biggest differences between SIRD 
and PSLQ. Furthermore, the matrix reducing method in 
SIRD is generalized Hermite reduction, which avoids LLL- 
type reduction. This is a difference not only between SIRD 
and HJLS, but also between SIRD and PSLQ because that 
(modified) Hermite reduction is not suitable to detect SIRs 
any more. And just the generalized Hermite reduction guar- 
antees the correctness of SIRD. 

4. PERFORMANCE RESULTS 

In theory, the costs of SIRD and the HJLS simultaneous 
relations algorithm (see [171 section 5]) are the same as in 
Corollary 13.61 in the worst case, whereas in practice SIRD 
usually needs fewer iterations. For vi — (11,27,31)"^ and 
V2 = (1,2,3)^, HJLS outputs (19,-2,-5)^ after 5 itera- 
tions while SIRD outputs (—19, 2, 5)"^ after only 2 iterations. 



No. 


n 


itrnjLS 


itrsiRD 


tujLS 


tsiRD 


1 


4 


15 


12 


0.047 


0. 


2 


4 


13 


9 


0.171 


0.016 


3 


4 


21 


19 


0.062 


0.015 


4 


5 


25 


20 


0.110 


0.016 


5 


5 


27 


43 


0.125 


0.016 


6 


5 


21 


14 


0.110 


0.032 


7 


30 


51 


21 


1.703 


0.422 


8 


54 


34 


9 


5.625 


1.265 


9 


79 


34 


40 


14.157 


4.422 


10 


97 


37 


5 


23.860 


5.375 


11 


128 


45 


6 


49.657 


11.141 


12 


149 


29 


14 


76.797 


18.063 


13 


173 


26 


2 


114.140 


25.000 


14 


192 


29 


2 


153.078 


33.641 


15 


278 


28 


8 


440.781 


102.860 


16 


290 


35 


6 


500.562 


118.578 


17 


293 


23 


7 


512.796 


123.265 


18 


305 


22 


4 


581.844 


137.672 


19 


316 


19 


3 


649.032 


147.796 


20 


325 


18 


2 


716.094 


159.813 



Table 1: Comparison of performance results for 
HJLS and SIRD 

Both the SIRD algorithm and the HJLS simultaneous re- 
lations algorithm when t = 2, i.e. detecting an SIR for two 
vectors, were implemented in Maple 13 by the first author. 
The tests were run on AMD Athlon'^" 7750 processor (2.70 
GHz) with 2GB main memory. 

The purpose of the trials in Table [T] is to compare the 
performances of HJLS and SIRD. n in Table [T] gives the 
dimension of the relation vector. itrnjLS and itrsiRD are 
the numbers of iterations of HJLS and SIRD respectively. 
The columns headed tHjLS and tsiRD give the CPU run 
time respectively of the two algorithms in seconds. 

The 20 trials in Table [T] were constructed by Maple's 
pseudo random number generator. The first 6 trials are for 
low dimension, and others for higher dimension. The results 
show that SIRD appears to be more effective than HJLS. 
In 18 out of 20 trials, the number of iterations of SIRD is 
less than that of HJLS. It is still true that SIRD usually 
needs fewer iterations than HJLS for more tests. This leads 
that the running time of SIRD is much less than HJLS. 
With n increasing, the difference between the efficiency of 
SIRD and HJLS is increasingly notable. On average, the 
SIRD running time is about 26.7% of the running time of 
HJLS. All these results are obtained under the condition 
that 7 = 2/V3+ 10"^". 

The Maple implementation and more tests are available 

from |http : //cid- 5ilbbl6a211c63a9b ■ skydrive . live . com/self . aspx/ . Public/sird . rar ] 

5. AN APPLICATION 

Any SIR detecting algorithm intervenes in many fields of 
application, such as Diophantine approximating, numerical 
constants relations finding, etc. In this section, we discuss 
how to find the minimal polynomial of a complex algebraic 
number from its an approximation by using SIRD. 

5.1 The MPF Algorithm 

We say that a complex number a is an algebraic number 
if a is a root of a non-zero polynomial in one variable with 
integer coefficients. The minimal polynomial of a is the 



unique primitive polynomial p{x) G Zlx] of least degree such 
that p{a) = 0. The degree and height of a are the degree 
and height of its minimal polynomial p{x) respectively. 

In this section, let a = a -I- 6/ £ C be an algebraic number 
with degree at most n, height at most H, where I — y/—T. 
Suppose we are given an approximation a to q such that 



max la — a 

l<i<n 



< e. 



(5.1) 



Is it possible to infer the minimal polynomial from the ap- 
proximation? Computer algebra system Maple has an LLL- 
based procedure, PolynomialTools : -MinimalPolynomial ( ) , 
for finding the minimal polynomial of an algebraic number 
from its an approximation, whose basic idea is from [271 1191 
120] . Applying SIRD, we shall give another affirmative an- 
swer, the following MPF algorithm, to the question above. 

Algorithm 4 (The MPF Algorithm). 



Input: an approximation q to q satisfying (|5.ip . a degree 

bound n, and a height bound H, e satisfying (|5.3p 
Output: the minimal polynomial of a. 
while 2 < i < n do 
v:=(l,a,--- ,aY 

Call SIRD with j — 2 producing an integer relation 
Pi = (Po,Pi,--- ,Pi)^ for V 



Pi :=the primitive part of "^'^^qPjx-' 

if height(pi) > 2"-^y/^irTTH then 

i := i + 1; goto Step[T] 
else return pi 
end while 



Remark 5. At Step [3] of MPF, pi is an SIR for vi = 
(l,Re(Q),--- ,Re{a')f and V2 = (0,Im(Q),--- ,Im(a"))^ 
when Im(a) 7^ 0. 

5.2 Error ControUing 

The main idea of our minimal polynomial finding (MPF) 
algorithm to determine the minimal polynomial of an alge- 
braic number from its an approximation is as follows: We try 
the value of i = 2, ■ • ■ , n in order. With i fixed, we call SIRD 
for detecting an exact integer relation pi = {po,pi, ■ ■ ■ ,Pi)'^ 
for V — (l,a,--- ,a')'^. Then pi(a;) = Yl^^oPi-'''' satisfies 
Pi{a) = 0, however, from which we can not decide whether 
Pi{a) is or not. Hence the most important problem is how 
to choose an appropriate e in (|5.ip such that Pi{a) — im- 
plies Pi{a) = 0. Before describing it in detail, we consider 
the following example. 

Example 1. Let a = 2 + \/3I. We know that the minimal 
polynomial of a in I.[x] is 7-4x + x'^. Let a = 2.000 + 1.732/ 
be the approximation to a with four significant digits. Hence 
(l.,2.,l.)^, V2 = (0., 1.732, 6.928)^. Feeding SIRD 

V2 after 2 



Vl 

vi, V2 as its input vectors gives an SIR for vi 
iterations. The corresponding matrices B are 



/ 2 



1 



-1 -1 



V 











/ 7 2 \ 







-4 -1 


1 ) 




^110/ 



It is obvious that the first column of the latter one is an 
SIR for Vl and V2, and corresponds to the coefficients of 
the minimal polynomial of a. However, if we take only 3 



significant digits for tlie same data, after 3 iterations SIRD 
outputs (1213, -693, 173)^, wfiich is an SIR for (l.,2., 1.)^ 
and (0., 1.73, 6.93)^, but does not correspond to the coeffi- 
cients of tiie minimal polynomial of a. For this reason, we 
have to appropriately control the error such that the output 
of MPF is correct. 

Lemma 5.1. Let f be a polynomial in 'L\x\ of degree n. If 
max;i<i<„ |a* — a* j < e, then |/(a) — /(q)| < e ■ n ■ height{f). 

Definition 5.2 (Mahler measure). For any polyno- 
mial g = X^i^oS*^' ^ ^[-^l '^f degree m with the complex 
roots zi, Z2, ■ ■ ■ , Zm we define the Mahler measure M(g) by 

m 

M{g) = \grr.\Y[ma.^{l,\zj\}. 

j=i 

The Mahler measure of an algebraic number a is defined to 
be the measure of its minimal polynomial. 

Lemma 5.3. (see '22, Lemma 3]) Let Qi,--- , be al- 
gebraic numbers of exact degree of di, - ■ ■ ,dq respectively. 
Define D — [Q(ai, • • • , Oq) : Q]. Let P G Z[xi, • • • , Xq\ have 
degree at most Nh m xh ft < h < q). If P{ai, ■ ■ ■ , Uq) 7^ 0, 
then 

\P{a,,...,aq)\ > ||P||i-^nM(«,)-^^^/'^\ 

h = l 

where AI{a) is the Mahler measure of a. 

This lemma gives a lower bound on \ P{ai, . . . ,aq)\ if P(cii, 
• • • 7^ for an arbitrary multivariate polynomial P £ 
Z[xi,--- ,Xq]. If we apply it to g{x) = X^I^oS'^' 1\x\, 
then we have 

Corollary 5.4. Let a be an algebraic number with ex- 
act degree no and g{x) = '^^gQiX^ £ Z[x]. Suppose both 
height{g) and height{a) are < H . If g{a) 7^ 0, then 

\g(a)\ > {m+l)-^"°-^^ ■ (no + 1)'"^ • Ji-t^+^o-i), 

where height{a) is the height of a 's minimal polynomial. 

Proof. For f{x) G Z[x] with degree n, we have Lan- 
dau's inequality: M{f) < \\f\\2 (e.g. see [HI p. 154]), 
height(/) < ll/lli < (n + l)height(/), and height(/) < 
II /II 2 < + lheight(/). This corollary easily follows from 
Lemma 15.31 and the three facts above. □ 

Next we investigate how to choose e to enable MPF to 
correctly return the minimal polynomial of a from d. We 
denote the exact degree of a by no(< n). For 2 < i < n, 
Step O in MPF gives a polynomial pi £ Z[x] with degree < i 
such that Pi{a) = 0. From Corollarv 15.41 we know that if 
Pi{a) / 0, then 

\p^(a)\ >{i + l)(^-"«' ■ (no + 1)-* • > M, (5.2) 

where M = (n + l)-(i"-i) ■ 7y-(2"-i), 

Theorem 5.5. Let a, 6t and M be as above, andp a poly- 
nomial in Z[x] with degree < n and height < H . Then there 
exist some e such that |p(a)| = implies p{a) = 0. 

Proof. Set deg(p) — i{< n). From Lemma |5. II we have 
|p(a)| = \p{a)-p{a)\ < e ■ i ■ H < e ■ n ■ H . Thus ife< 
then |p(a)| < M. From ((5^ it follows that p(a) = 0. □ 



If we substitute 2'^-'^y/nTTH for we have 

Corollary 5.6. Let a and a be as above and 

e < 2-''"'+^"(n + l)~5"Ji-2". (5.3) 

Then for i from 1 ton, an integer relation for {1, a, ■■ ■ ,a^)'^ 
with height < 2""^i/n -I- IH is also for (1, q, • ■ • , a')^. 

5.3 Correctness and Cost of MPF 

Assume that the degree of a is no and that e satisfies (|5.3p . 
When 2 < j < no, there exists no relation for (1, a, ■ • ■ , a'), 
which, combined with Corollarv 15.61 means that Pi{x) must 
satisfy the condition in StepUJof MPF and then go into next 
iteration. When i = no (< n), we know that the coefficients 
of the minimal polynomial of a form an integer relation for 
(l,Qf, •■• ,a"°), whose height < H, hence Euclidean norm 
< y/no -\- IH. This implies that (1, a, ■ • ■ , a"°) has also an 
integer relation with Euclidean norm < \/no 4- IH . From 
Theorem 13.21 we know that the height of the relation SIRD 
detected will < 2"~^^n -I- IH. Thus the relation detected 
by SIRD when i = no will never satisfy the condition in 
Step [4] and corresponds an integral multiple of the minimal 
polynomial of a. Hence the correctness of MPF follows. 

From 1)5. 3[l we have log e G 0{n^ -\- n log H). Thus we can 
give another answer to Blum's and Zhang's question without 
using LLL lattice reduction algorithm. 

Theorem 5.7. Let a be an algebraic number and let n 
and H be upper bounds of the degree and height of a respec- 
tively. Suppose we are given an approximation a to a such 
that maxi<i<„ |a* — ti'j < e. Then the minimal polynomial 
of a can be determined m 0{n^ -f- n** log i?)) arithmetic op- 
erations on floating-point numbers having 0{n^ + nlog_H")) 
bit- complexity. 





Digits 


Complexity 


KLL^ 


0{n^ + n\ogH) 


0{n^ -i-n" log H) 


Just 18 


0{n^ + n^ log H) 


0{n'^ log n + n** logH) 


QFCZ 25 


0{n^ + n\ogH) 




MPF 


0{n-' + n\ogH)) 


0{n''-\-n'' log//)) 



Table 2: Comparison of different minimal polyno- 
mial finding algorithms 



Table [2] gives a comparison of the digits and complexity 
of 4 different minimal polynomial finding algorithms in the 
worst case. Since the algorithm in 25 can only find the 
minimal polynomial of a real algebraic number, we don't 
compare the complexity with it. It seems that a lower com- 
plexity can be achieved by using some new type LLL algo- 
rithms, such as L^ [24] and H-LLL [23], but when we apply 
these new algorithms to find the minimal polynomial we 
have to choose e as in a similar formula with 1)5. 3[) . Thus 
multiple precision arithmetic is inevitable. 

Example[l\ (con.). For a = 2 -\- y/3I, its minimal polyno- 
mial 7~'ix-\~x'^. Set n = 2 and H = 7. Computing the error 
tolerance as in equation (|5.3|l gives e < 583443"^. Corollary 
l5.6l imDlies that [— logj^Q 583443"^ J = 5 correct decimal dig- 
its are sufficient to guarantee the output is correct. This 
example also illustrates that e in (|5.3p is only a sufficient 
condition on error controlling, but not a necessary one. 



6. CONCLUSION 

The number of iterations and the cost of SIRD algorithm 
are related to the parameter 7. For Vi — (86, 6, 8, 673)"^ and 
V2 = (83, 5, 87, 91)^, if we choose 7 = 1.16 then SIRD out- 
puts (—215, 402, 159, 22)"^ after 12 iterations, however, if we 
choose 7 = 5, SIRD outputs (93, 364, 93, -14)^ after only 6 
iterations. In future work we expect to find the best choice 
for 7. Additionally, how to choose the digits such that SIRD 
under floating-point arithmetic flnds an exact SIR is also in 
our interests. Finally, we see that the MPF algorithm can 
be used to factor / in Z[x] like this: Solve an approximation 
root with accuracy satisfying equation (|5.3p . and call MPF 
for finding its minimal polynomial which corresponds an ir- 
reducible factor of /, and then repeat the two steps until 
/ has been factored completely. It is symbolic- numeric and 
different from traditional algorithms based on Hensel lifting. 
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